%%Jicamarca ISR data From Patrick
%read in, calculate average drift 
data = load('c:\manoj\projects\longp\jicamarca_isr.dat');
fday = data(:,1)/(24*3600)+datenum(1970,1,1); %CEDARWEB time stamp is in seconds from 1.1.1970
K = unique(data(:,1));
drift = zeros(size(K));
for i = 1:length(K),
L = data(:,1) == K(i);
mydata = data(L,[6,9]);
L = mydata(:,1) <400; % The drift data are for different altitudes. here selected all the altitudes less than 400 km
drift1(i) = nanmean(mydata(L,2));
L = mydata(:,1) <200; % The drift data are for different altitudes. here selected all the altitudes less than 200 km
drift2(i) = nanmean(mydata(L,2));
end;

%save c:\manoj\projects\longp\jicamarca_height_averaged drift1 drift2 FF;
%%
FF = unique(fday);

L = isnan(drift2);

drift = drift2(~L);
FF(L) = [];

% drift = drift2(~L).*24.366*1e-3; %mV/m
drift = drift2(~L).*24.366*1e-3; %mV/m

frq= fftfrq(256,1/(15*60));%frequenxy in Hz

prd = 1./(frq*24*3600); % period in days


for i = 1:length(frq),
    c(i,:) = inversion(FF'-datenum(2000,1,1),drift,prd(i));
end;

nfrq=1;
for i = 2:-0.01:0.02,
    c(nfrq,:) = inversion(FF'-datenum(2000,1,1),drift,i);
    nfrq=nfrq+1;
end;
